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Abstract 


This article deals with an efficient approximation method named succes- 
sive complementary expansion method (SCEM) for solving singularly per- 
turbed differential-difference equations with mixed shifts. It is compared 
with the method of matched asymptotic expansion (MMAE) and the pa- 
rameter uniform upwind finite difference scheme for solving such a model. 
The comparison shows, unlike the MMAE, the SCEM method requires no 
matching procedure. It requires less computation when compared to the 
upwind finite difference scheme on the Shishkin mesh. The error analy- 
sis is carried out to prove the robustness of the method. Some numerical 
experiments are provided, which show the effectiveness of the proposed 
method. 
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1 Introduction 


It is noteworthy to observe the wide use of differential equations involving 
small parameters in a variety of disciplines from physics, chemistry, eco- 
nomics to biology. A differential equation, in general, is said to be singularly 
perturbed if its higher-order derivative gets multiplied with a small parameter 
é (0 <é <1) known as the perturbation parameter. In many mathematical 
models, these types of equations contain some small delay and small advance 
terms, which make them singularly perturbed differential-difference equa- 
tions (SPDDE). Due to the presence of the perturbation parameter, these 
problems exhibit boundary or interior layers. The development in solutions 
to these problems both analytically and numerically is not a settled direc- 
tion in mathematics. The massive use of these equations in many fields of 
science and engineering leads many researchers to solve such problems more 
efficiently. These equations usually arise in control theory [12], generation 
of action potential in nerve cells by random synaptic inputs in dendrites, 
modeling of activation of a neuron, and many more [25, 27]. 

Over the last few decades, many researches have been done on singularly 
perturbed differential equations. One may refer to [1, 6, 5, 13, 14, 24] and 
the references therein for different kinds of numerical methods. To mention 
a few, Kadalbajoo and Sharma [10, 9] surveyed various asymptotic and nu- 
merical methods for the solution of such problems. Lange and Miura [15, 16] 
gave an asymptotic approach to solve SPDDE along with the turning point 
problems. Sirisha, Phaneendra, and Reddy [26] used the idea of domain 
decomposition and introduced the mixed finite difference method. Using 
post-processing and grid distribution techniques, Mohapatra and Natesan 
(19, 20, 21] formulated uniformly convergent numerical methods. Duressa 
and Reddy [8] constructed a domain decomposition method and used a ter- 
minal boundary point to the domain for decomposing it into two regions. 
Melesse, Tiruneh, and Derese [18] used an initial value method for solving 
these problems. A collocation method was proposed to solve these bound- 
ary value problems using a modified B-spline basis function by Arora and 
Kaur [2]. Rao and Chakravarthy [23] used an exponentially fitted tridiago- 
nal finite difference method. Recently, Cengizci, Natesan, and Atay [4] used 
an asymptotic numerical hybrid method for the singularly perturbed system 
of two-point reaction-diffusion boundary-value problems. Mushahary, Sahu, 
and Mohapatra [22] used a finite difference scheme for solving SPDDE. 

Our main objective in this paper is to prove a successive complemen- 
tary expansion method (SCEM) as an efficient alternative to the method of 
matched asymptotic expansion (MMAE) for solving SPDDE. This method, 
in general, applies to all the singularly perturbed problems that can be ap- 
proximated by the MMAE. The main idea of SCEM is to use a correction 
term for the approximation of the inner layer to complement the solution for 
the entire domain. The quality of the approximation is improved by itera- 
tively using corrections to the new terms for the inner layer. The ability to 
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provide a uniformly valid approximation without any matching procedure is 
the main advantage of using SCEM over MMAE. 


In this work, some SPDDEs are taken into consideration, which exhibit 
boundary layer behavior at the left-end of the interval. These problems 
contain mixed shifts, that is, both the delay and the advance terms. The 
complete description of the problem is given in section 2. At first, the model 
is approximated asymptotically by the MMAE, and then we have used an 
efficient asymptotic method as an alternative of widely used MMAE, that 
is, the SCEM. This method was first introduced by Cousteix and Mauss 
in [7, 17]. We have also approximated these problems by an upwind finite 
difference scheme on the Shishkin mesh. The complete working principle of 
the parameter uniform numerical scheme, MMAF, and SCEM, along with the 
use of SCEM over MMAE, has been described in section 3. In the support 
of the above-said methods, some numerical examples with their results are 
discussed in section 4. The complete discussion with the conclusion is given 
in the final section. 


2 Description of the problem 


Consider the following SPDDE with mixed shifts terms: 


+e(x)y(x) + d(z)y(z +) = f(z), 0<2<1, 
; on —d<2<0, (1) 


; on l<a<1+4+y, 


where a(x), b(x), c(x), d(x), f(x),o(x),y(@) are bounded and continuously dif- 
ferentiable functions on (0,1). Here, ¢ is the singular perturbation parameter 
with 0 <<e <1. The parameters 6 and 7 are the delay and advance terms, 
respectively, with 0 < 6 = o(€);0 < 7 = o(e). Now for the terms of (1) con- 
taining the delay and the advance operators, we shall use the Taylor’s series 
approximation in the neighborhood of the point z, that is, 


y(a — 6) © y(x) — dy'(a), 


re ; (2) 
y(x +n) © y(x) + ny’ (2). 


The application of (2) converts (1) to the following form: 


fe + a(x)y"(x) + B(a)y(@) = fz), O<a<1, 
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where a(x) = a(x) — 6b(a) + nd(x) and B(x) = W(x) + c(x) + d(x). Assume 
a(x) > A > 0 and G(x) > 0 for some constant A. Since the parameters 6 and 
7 are of o(€), so the use of Taylor’s series approximation and the conversion 
of (1) to (8) is valid. 


3 Description of the methods 


3.1 Asymptotic approximation 


Any two functions w(z,¢) and €(#,¢) are said to be asymptotically identical 
to p(e) if their difference is less than p(¢) in the domain Q, 


W(w,€) — €(@,€) = o(p(e)). (4) 


Here, p(e) is an order function and €(,¢) is the asymptotic approximation 
of w(x,¢). The general form of asymptotic expansion is 


b(a,e) = S_ di(e)yi(a,), (5) 
i=0 
where 6;(€), i =0,1,..., is an asymptotic sequence and w;(a,¢) is of o(p(e)). 


If d;(€) satisfies d;41(€)= o(d;(e)), then the expansion is called generalized 
asymptotic expansion. If the expansion is of the form 


(a, €) = Bow = >> (4! (a, 2), (6) 


1=0 


then it is said to be regular asymptotic expansion. Here, Eo is an outer 
expansion operator. So, we conclude w(a,<¢) — Eow(#,¢) = o(p(e)). 


3.2 Method of matched asymptotic expansion (MMAE) 


In many cases, the functions are not regular in the domain of consideration 2. 
So, (6) is only valid upon one finite part of the region 0. Let that restricted 
region be named N59 C 2. Now, the other part of the region (inner region), 
that is, 2, is located near the boundary layer where the solution undergoes 
a rapid change to satisfy the boundary conditions. For the inner region, 
we need to introduce a new variable(stretching variable or boundary layer 
— 20 
(a 


solution occurs and €(¢) is the order of the thickness of the boundary layer. 


variable) t=" Here, zo is the point where the rapid change of the 
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Now, similarly as in (6), we can construct an asymptotic expansion for the 
function in Q,, 


b(a,e) = Erb = 98. (ew (a, e). (7) 


Both the expansion operators Ey) and Fy have the same order, that is, o(p(e)). 
So, we have w(a,¢) — E,W(%,¢) = o(p(e)). The valid approximation of the 
solution thus can be written as 


Wmmae(X, €) = Eow(x,€) + E,y(@, €) ~~ Eo E(x, €). (8) 


In MMAE, we calculate two different results separately with two different 
variables for both the inner and outer layers. Then, to obtain the uniformly 
valid approximation for the entire domain, we need to subtract the region 
where both the approximations overlap. For this purpose, we need to develop 
a matching procedure by the limiting process [3]. 


3.3 Successive complementary expansion method 
(SCEM) 


As the functions that we consider are not regular in the entire domain Q, after 
calculating the solution for the outer region Qo through the expansion Eo as 
in (6), we need to extend the solution through a complementary function for 
the entire domain. This complementary function depends upon the variable 
z. The general form of the asymptotic expansion in SCEM can be formulated 
as 


Wscem(x, £; €) = Dae we" (x, €) + Yi", €)). (9) 


Here, ~?"'(x,€) is the same as the outer region solution that has been ob- 
tained by MMAE before. The outer solution depends upon x but not upon 
€ but the complementary function depends both upon x and «. As in the 
MMAE method, SCEM does not require any matching procedure. 


3.4 Finite difference scheme 
In this section, we first construct a piecewise uniform Shishkin mesh. The 


Shishkin mesh is distinguished from other piecewise uniform meshes in terms 


of its choice of the transition parameter 7 defined as 7 = min (3. -lnN ), 
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Now, we divide each of the intervals [0,7] and [7,1] into N/2 equal sub- 
intervals. The theory of Shishkin mesh demands as N — oo the number of 
sub-intervals in each interval should be bounded below. 


_ jet), 15 .N/2, 
ao" = oe =f, 45, oo 


2(1 — 7) 


The coarse part of Shishkin mesh has spacing H = and the fine 


2T 
part has spacing h = —. Now, let us denote h; = x; — x;_,. The forward 
and central difference operators are denoted as 


Yin -—V% 


2 Y; —¥Y; Y; -— Yi-1 
DtD 41 
hig 


Dt= = 
hit hisa higt hj 


Now (3) takes the form 


eD* D-Y, + a(x)D*Y; - B(2)¥; = f(e), 
y0)=YO), yQ@=YQ), (11) 


for any mesh function Y(a;) = Y;. So, for the domain, in (10), we get the 
following system of equations: 


P,Y;-1 See a =P, (12) 
where 
2€ 
PO SS 
: ha(he + higi)’ 
2e 2e a(a) 
Po = + B(x), 13 
hi(hi t+ higt)  higi(hit hig) igi Bm) (13) 


P= 2€ x a(x) 
higi (hi t+ higt) igi 


The matrix associated here is an M-matrix. For detailed analysis of the 
proposed scheme one can refer [11]. 


4 Results and discussion 


Even though much work has been done by the MMAE method, there is no 
specific theory that can describe the method in a proper way. So, we need to 
discuss the working principle of both MMAE and SCEM by some illustrative 
examples. The proposed SCEM is found efficient in solving examples of the 
type (1), which possess the left-end boundary layer. The approximated exact 


IJNAO, Vol. 12, No. 1, (2022), pp 55-72 


Asymptotic and numerical methods for solving singularly perturbed ... 61 


solution of SPDDE of the form (1) with constant coefficients is given by 


y(x) — Ce + Coe (14) 
where 
1 —eM2 
C= eM1 — @M2? C2=1-Ci, 
—a(x) + /a(x)? — 4B(ax —a(x) — ,/a(x)? — 4B(x 
my = ate) + VP TBE) y, _ ale) — Vora), 


Example 1. Consider the following problem with a delay term: 


{ ey! (x) + y'(x) + 2y(x — 6) — 3y(x) = 0, 
yj)=1, -d<a2<0, y(l)=1. 


The solution of this problem possesses a rapid change near the point x = 0 
when ¢ + 0+. This region is said to be the boundary layer or the inner layer. 
So, this possesses a left-end boundary layer. The other region, which is far 
away from the boundary layer, is said to be the outer region. In case of 
the outer region, there is no record of any unusual behavior of the solution. 
We will solve this problem employing MMAE and SCEM and the proposed 
numerical scheme. For the given example using Taylor’s series expansion of 
the form (2) in the neighborhood of x, we have 


* 1 * E 
e*y" (x) + y/(x) — Toa ¥) = 0, where &* = Toa (15) 


Now, the equation has converted to a second order singularly perturbed prob- 
lem, which can be solved by both MMAE and SCEM. 


4.1 For MMAE 


In MMAF, for the outer region solution, which is far away from « = 0, we 
take « = 1. Now, the asymptotic approximation for the outer region can be 
written as 


y(x) © yo(x) + €*yr(x) + (€*)*yo(x) ++. (16) 
Substituting (16) in (15), we reach at 
e*[yg (x) + e* yi (2) + (€*)?yo(x) +--+] + [yo(a) +e" (x) + (17) 


(*)usla) ++] = (755) ule) + evn le) + Punto) ++] =0 


For ¢ = 0, model (17) is reduced to 
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solo) ~ (7255 }uola) =0. (18) 


1— 26 
A € R and imposing the outer boundary layer condition (at 2 = 1) we get 
our required outer solution as 


yo(x) = exp (5). (19) 


To obtain the solution for the inner layer, we need to introduce the stretching 
variable = =. The solution obtained in the inner region that depends upon 
the variable Z, is denoted by Y(%). By the use of this stretching variable, we 
are able to stretch the thin layer near the boundary. Now by applying the 
chain rule, we get 


The solution of (18) is easily found to be yo(a)=A exp ( = ), Here, 


ai ea 2 
dx edz 20) 
Using (20) in (15), we have 
PY) dY (&) 1 
* 1 *\—1 Ny 
(ey ar + eye a apY@=0. Qt) 
Now, multiplying (21) with «*, we have 
a@Y(z)  dY(z) e* 
Y(z) =0. 22 
dz? ic Co seal 


Equation (22) is a regularly perturbed differential equation, and to obtain the 
solution for the inner region, we need to use the asymptotic approximation 
of the form 


Y(z) © ¥Yo(z) + e*Vi(z) + (e*)?¥o(z) +-->. (23) 


We are only doing the expansion for the first term. So, we can make our 
calculations for ¢ = 0. Now, we have 


Yo (£) + Yo(#) = 0. (24) 


On solving (24) and using the boundary condition for the inner layer (at 
x = 0), we arrive at 


Yo(Z) =1+C(e* — 1). (25) 


Clearly, there are no other conditions, so that we can obtain the value of the 
constant C’. Here is the situation where the matching procedure of MMAE 
works. We have already found out two different solutions, that is, (19) and 
(25) for two different regions, but, it is obvious that they belong to the same 
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approximation. Followed by this matching idea [3], we conclude 


lim yo(2) = lim Yo(z). (26) 
x—>0 @—-0o 


1 
Thus, we get C = 1 — exp (ss): So, the solution of our inner region 


becomes 


Vota) = exp ( =) ex (559 z) +e (55-7). (27) 


As we have calculated the solutions for both the outer and inner regions, now 
we can find out the composite solution to our equation. For this, we need to 
add the solutions of both the regions and subtract the common limit between 
them. As per (8), we have 


y © yo(x) + Yo(Z) — yo(0*), (28) 


Using (28), we obtain the composite MMAE solution as 


S x—1 a —2x 1 x (29) 
Ymmae Ce exp 1 — 26 exp € exp 26 _ i e * 


4.2 For SCEM 


The uniformly valid composite expansion in general form of SCEM can be 
written as 


Yocen(a;%,£) = > _ d:(e)[ys(a) + ¥i(2))- (30) 
i=0 
As we are only interested in solving our problem for x = 0, so we get 
Yoscem(a, ©, €) = [yo(a) + Yo(@)]. (31) 
Here, yo() is the same outer region solution that has been obtained for 
MMAE in (20). The outer region solution remains the same as the solution 


here showing no unusual behavior in accordance with the boundary condi- 
tions. Using (20) in (31), we get 


—1 
Yscem = €XP (5) + Yo(Z). (32) 


Substituting value of Yscem in (15), we have 
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1 
2 Vossen (2,2, 9) + isan (z, a, é) = T— 26 Yoscem (tit €) = 0. (33) 
This follows 
es z—1 7 d z—1 
E [ex (5) +n@] +2] ov (ES) (34) 
1 x—1 
Yo (& — Yo(z)| = 0. 
o(0)] - [ew (5) + ¥0@] =0 
Using the chain rule for the solution of Yo(Z) in terms of x, we have 
a c—1 1 1 
* Se pa oe aay = pee 
“|(=a) © (a) + eee O|+ (Tw) © 
x—1 1 1 x—-1 
ha) a) vate) |. 0 
=e (; = =i) +S tc) ay, exp ( = =i) " o(@| 
Multiplying e* in (35), the resulting equation will be 
42 1 x—1 L¥i"(z) +e" 1 ag z—1 
= =95| ieee) Re PS | tags aay 
+(e) —et|—_]] exp (271) +m] =0, (38) 
— em liga yl a 


which is regularly perturbed ordinary differential equation, and for e* = 0, 
it will be reduced to 


Yo (€) + Yo(®) =0. (37) 
Now, the solution of the above equation is given by 
Yo (Z) _ Cy + Cye*. (38) 


For a uniformly valid SCEM solution, Yoscem must satisfy the following 
boundary conditions: 


t= 0 > E = = =0> Youcem(0,0,£) = yo(0,2) + ¥o(0,0,2)=1 (39) 
=> Yo(0,0,¢) = 1 — yo(0,¢), 


8 

a 

al 
| 


1 1 
L Yoscem(1, 276) = yo(1, €) + Yo(1, 276) =1 
1 
= Yo(1, —,€) = 1 —yo(1,¢). 


Using the boundary conditions in (38), we get 


IJNAO, Vol. 12, No. 1, (2022), pp 55-72 


Asymptotic and numerical methods for solving singularly perturbed ... 65 


(c) « = 0.009 (d) e = 0.0004 


Figure 1: Comparison of solutions with different values of ¢ for Example 1. 


= 7 i 
= 1 5 Ci = —C2 exp (=). (40) 
exp (=) —1 ' 
e* 
Now, using (40) and (39) in (32), we get the composite SCEM solution as 
x—-1 exp(7=3;) — 1 —1 
scem = €X ex ae 
GY ap )—-1 | Se 
exp(7—35) -—1 —2£ 
+ ex : Al 
aba [or (S = 
After getting the approximation to the exact solution by both SCEM and 
MMABE, it is quite obvious that SCEM is more adaptable. The immense 
advantage of SCEM is that it gives the uniformly valid approximation to the 
exact solution without any matching procedure. The error in comparison 
with the approximated exact solution for both SCEM and MMAE is given 
in Table 1 for different values of ¢. The result is proved to be valid for any 


value of N. Here for the computational purposes, the nodal points and the 
values of 7 and 6 are kept fixed to N = 32 and 7 = 0.le = 0. 


Example 2. Consider the following problem having an advance term: 


ae + y' (a) — 3y(x) + 2y(x +) =0, 
yO=1, y(x)=1l<ax<1l+n. 


This problem has a rapid change near z = 0. Now using Taylor’s Series 
expansion in the neighborhood of x, the above problem reduces to 


e*y" (x) + y'(z) (x) = 0, where e* = Z 


~ 1429" 14+2n 


(42) 
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E L®errorin MMAE | L®errorinSCEM | Error in Upwind 
0.0001 0.00003676 0.00003676 0.04450943 
0.0005 0.00018374 0.00018374 0.04428418 
0.0010 0.00036728 0.00036728 0.04401059 
0.0050 0.00182533 0.00182533 0.04209958 
0.0100 0.00361593 0.00361593 0.04025378 
0.0500 0.01598307 0.01598307 0.03472166 
0.1000 0.02552850 0.02556347 0.03481582 
0.3000 0.02853413 0.03217078 0.03427254 
0.4000 0.06644741 0.03539967 0.03541439 
0.6000 0.15664433 0.04561360 0.03723168 
0.8000 0.24352992 0.05726515 0.03825982 
1.0000 0.32059406 0.06915099 0.03877738 


a “| | eee Pg 
\ 3 we ‘ \ we 
(a)e=0.1 (b) € = 0.05 
a =e 
‘seen La a 
(c) ¢ = 0.006 (d) ¢ = 0.0009 


Figure 2: Comparison of solutions with different values of ¢ for Example 2. 


The one-term SCEM and MMAE approximations for (42) are given as follows: 


1 —1 —2x£ i: x—1 
mmae = |1— exp {| ——— ] | ex ex ; 
v se a + 2n . e* el + 2n 


seg, a [OPTED =") oy (2) 4 [PPD] oop (2 
ia exp(=) -1 es exp(=})-1 e* 


e* 


Example 3. Consider the following problem having mixed shifts: 
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E L®errorin MMAE | L°errorinSCEM | Error in Upwind 
0.0001 0.00003676 0.00003676 0.04448748 
0.0005 0.00018367 0.00018367 0.04417452 
0.0010 0.00036699 0.00036699 0.04379147 
0.0050 0.00181786 0.00181786 0.04101098 
0.0100 0.00358632 0.00358632 0.03808839 
0.0500 0.01533655 0.01533655 0.02519028 
0.1000 0.02377291 0.02379587 0.01688285 
0.3000 0.01783685 0.02658575 0.01728136 
0.4000 0.04058107 0.02660582 0.02334674 
0.6000 0.09131635 0.02829255 0.03160623 
0.8000 0.13551419 0.03009785 0.03671486 
1.0000 0.17029574 0.03133850 0.04028387 


Table 2: Numerical results with different values of ¢ for Example 2. 


This problem has a rapid change near x = 0. By the use of Taylor’s Series 
expansion for both the delay and shift terms as in (2), we have 


6 


* * E 
cy'(a)+¥(@)- | ares |ule)=0, where e* = 


— 1+264+7 


The one-term SCEM and MMAE approximations are given by 


= |1 —ex =e ex on + ex Sie= 1). 


—6 —6 
exp 1 —1 exp 1 = 
ree | ( 1$25+q) | 2 ( ) | ( T2079) | é ( =) 


exp(=) —1 e exp(=)—1 ef 


6(a — 1) 
14+26+7]" 


+ exp | 


The comparison of errors of this example for both the proposed methods and 
the method prescribed in [10] for N = 8 and N = 256 are given in Tables 
4 and 5, respectively. These comparison results confirm the efficiency of the 
proposed method. 
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(d) ¢ = 0.0007 


Figure 3: Comparison of solutions with different values of e for Example 8. 


E L®errorin MMAE | L°errorinSCEM | Error in Upwind 
0.0001 0.00022006 0.00022006 0.05991807 
0.0005 0.00109602 0.00109602 0.05951004 
0.0010 0.00218145 0.00218145 0.05900561 
0.0050 0.01050670 0.01050670 0.05517982 
0.0100 0.02012147 0.02012147 0.05086176 
0.0500 0.07770075 0.07770075 0.03679054 
0.1000 0.12342561 0.12345914 0.03402622 
0.3000 0.16404633 0.18962499 0.01593674 
0.4000 0.16021248 0.19825547 0.01255964 
0.6000 0.15178975 0.20931461 0.00894834 
0.8000 0.21056755 0.21768491 0.00694484 
1.0000 0.26983417 0.22335443 0.00562079 


Table 3: Numerical results with different values of ¢ for Example 3. 


€ L*errorin MMAE | L®errorinSCEM | Error in [10] 
1-* 0.10596702 0.10597709 | 0.12011566 
10-2 0.01882187 0.01882187 | 0.18727108 
10-8 0.00209815 0.00209815 | 0.20429729 
10-4 0.00021228 0.00021228 0.20614146 
10-7 0.00002125 0.00002125 | 0.20632746 
10-® 0.00000212 0.00000212 | 0.20634608 


Table 4: Comparison of maximum error with N = 8 for Example 3. 
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€ L®errorin MMAE | L®errorinSCEM | Error in [10] 
10 * 0.10685019 0.10685019 | 0.00775036 
10-4 0.01976044 0.01976044 | 0.00799076 
10-3 0.00218116 0.00218116 | 0.00963304 
io-* 0.00022045 0.00022045 | 0.00984236 
10-° 0.00002206 0.00002206 | 0.00986365 
10~° 0.00000220 0.00000220 | 0.00986578 


Table 5: Comparison of maximum error with N = 256 for Example 3. 


5 Conclusion 


In this work, the well-known MMAE and a parameter uniform numerical 
scheme, that is, the upwind scheme in Shishkin mesh were compared with 
the proposed relatively new method, that is, SCEM, for solving singularly 
perturbed differential-difference equation with mixed shifts. It was observed 
that both MMAE and SCEM provide highly accurate approximations to the 
exact solution for comparatively small values ¢. Similarly, comparison results 
confirmed that the numerical scheme, which gives better results for higher 
values of V, at the same time, MMAE and SCEM give higher-order accuracy 
for comparatively smaller values of N. It is observed that we do not require 
any matching procedure in SCEM as required for MMAF, and the boundary 
conditions are satisfied exactly for all values of the perturbation parameters. 
These properties make SCEM more flexible and a better alternative while 
solving singularly perturbed differential equations. 
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